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Abstract 

We introduce a new stochastic model for the variations of asset 
prices at the tick-by-tick level in dimension 1 (for a single asset) and 2 
(for a pair of assets). The construction is based on marked point pro- 
cesses and relies on linear self and mutually exciting stochastic inten- 
sities as introduced by Hawkes. We associate a counting process with 
the positive and negative jumps of an asset price. By coupling suitably 
the stochastic intensities of upward and downward changes of prices for 
several assets simultaneously, we can reproduce microstmcture noise 
[i.e. strong microscopic mean reversion at the level of seconds to a few 
minutes) and the Epps effect (i.e. the decorrelation of the increments 
in microscopic scales) while preserving a standard Brownian diffusion 
behaviour on large scales. 

More effectively, we obtain analytical closed-form formulae for the 
mean signature plot and the correlation of two price increments that 
enable to track across scales the effect of the mean-reversion up to the 
diffusive limit of the model. We show that the theoretical results are 
consistent with empirical fits on futures Euro-Bund and Euro-Bobl in 
several situations. 
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1 Introduction 



With the availabihty of a huge number of quahty high frequency data, there 
is a fast growing hterature devoted to the modelhng of intra-daily asset prices 
behaviour. Since Bacheher and the seminal work of Black and Scholes, 
most popular models at a coarse time scale ~ for daily data, say - are 
Brownian diffusions, see for instance the classical textbooks of Musiela [25] or 
Bouchaud and Potters [S] and the references therein. In particular, diffusion 
models aim at describing more or less faithfully the volatility dynamics, 
characterized by stylized facts such as volatility clustering or leverage effect 
[8] . A key issue that naturally emerges when one studies high frequency data 
is the problem of improving volatility estimation or covariance estimation 
between two asset returns over a given time period, thanks to the massive 
amount of data available at such scales nowadays. The discrete nature of 
time trade arrivals and of price variations (prices are point processes living 
on a tick grid), the presence of so-called microstructure noise (described as 
strong mean reversion effects at small scales) makes this question highly non 
trivial. At a very high frequency, prices variations are also characterized by 
well documented stylized facts like the signature plot behaviour and the 
Epps Effect [Hj. 

1.1 High frequency volatility: microstructure noise 

If X{t) stands for the price of some asset at time t (defined indifferently as 
the last traded price or the mid-price between best bid and best offer in the 
order book), the signature plot can be defined from the quadratic variation 
of X{t) over a time periode [0, T] at a scale r > - the so-called realized 
volatility - as 

1 

C{r) = -Y,\X{{n+l)T)-X{nT)\'. (1) 

n=0 

The microstructure noise effect manifests through an increase of the ob- 
served daily variance when one goes from large to small scales i.e. in the 
limit T — 7- (see e.g. Fig. [5|^b)). This behaviour is different from what one 
would expect if the data were sampled from a Brownian diffusion, in which 
case the function plotted in Fig. El^b) should be flat. From the perspective 
of statistical estimation, this leads to a simple paradox. On one side, the 
smaller r, the larger is the dataset that can be use to estimate the volatility. 
However, how should one be using high-frequency data in order to obtain 
better estimates of the volatility, since the realized volatility ([T]) is not stable 
as r decreases. 

In the literature the most popular approaches attempts to model mi- 
crostructure noise with the concept of latent price. One starts with a Brow- 
nian diffusion X(t) defined as an efficient price, which is latent, in the sense 
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that it cannot be observed directly. Instead the practitioner has only access 
to a noisy version X{t) of X{t) that accounts for microstructure noise. The 
most successful version - as far as mathematical development is concerned 
- is the additive microstucture noise model, introduced in 2001 by Gloter 
and Jacod \15\ [T6] and in the context of financial data by Ait-Sahalia et al. 
[2t I35[ [M]. Given a sampling scale r, one rather observes 



where the microstructure noise term ^n,T satisfies E[^„^,-] = for obvious 
identifiability conditions. (Hereafter, E[-] denotes the expectation operator.) 
The goal is then to separate the noise from the true signal X{t), from which 
a classical volatility estimator can be performed. This has raised a vast 
research program over the last decade, mostly covered by econometricians 
and statisticians, see [IllMIIIlEIllMlEElllElElEOllSlEllESllig and 
the references therein. Whereas representation ([2]) produces an elegant pilot 
model to describe microstructure noise effects at the scale of a few minutes, 
it cannot faithfully reproduce the data as they are observed on a microscopic 
scale of a few seconds: for instance, the discreteness of price changes is left 
out and the mathematical artefact of forcing C(t) to explode when r — )• 
becomes unavoidable. 

1.2 High frequency correlations : the Epps effect 

Another important feature concerns the way different asset price movements 
are correlated. Stong correlations between asset returns at daily or larger 
time scales usually exist, and modern portfolio theory precisely relies on such 
dependences. One expects that coarse scale correlations originate from in- 
traday strongly correlated movements, and this naturally raises the question 
of how to get better estimates of such correlations by using high frequency 
data. Again, an analogous paradox as for volatility estimation under mi- 
crostructure noise arises. In Fig. [6l we display the last traded ask prices of 
futures Bund lOY and Bobl 5Y over a period of a few hours, prices incre- 
ments are clearly correlated in their intra-day variations. If Xi[t) and X2{t) 
are the prices of two assets, a correlation coefficient estimator over a time 
period [0, T] can be naturally defined from high frequency price increments 
as 




(2) 




where 



Ci2(r) = - ^ + l)r) - Xi(nr)] [X^iin + 1)t) - X^inr 



■)] 



n=0 
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and Ci(r) and C2(t) denote the realized volatility of Xi{t) and X2{t) re- 
spectively at scale r. In Fig. El^c), we plot p(r) as a function of r for the 
couple Bund/Bobl. It is striking to observe that the correlation coefficient 
is an increasing function of the time resolution and that correlation almost 
vanish at a very high frequency. This phenomenon, first reported by Epps 
[11], is the so-called Epps effect. Few approaches however address the Epps 
effect in the literature. 

Besides the problems related to the estimation of coarse scale asset prop- 
erties, we see that modelling high frequency price dynamics is a source a 
many challenging questions. Other approaches consist in defining "fine" 
scale models and addressing directly the price dynamics at the tick level, see 
for instance Engle and collaborators, |13| [T2] who introduce the ACM-ACD 
model or Bauwens and Hautsch [7], Bowsher [9] who construct intensity 
based point processes for modelling high frequency data. In this context, 
Hewlett [18j introduces a model for measuring trade arrivals that is formally 
very close to our construction. However, Hewlett's approach is mostly order 
book oriented: it focuses on the imbalance properties of an asset and leaves 
out the directional behaviour of the prices. In the same way, the aforemen- 
tioned literature that uses point processes is essentially specialized to high 
frequency data description, or focuses on trades arrivals dynamics. In par- 
ticular, the questions related to intermediate or asymptotic behaviour of the 
statistical price properties as the resolution scale varies are not considered. 
We plan to address in the present paper this next logical step. 

1.3 Objective and content of the paper 

We introduce a "fine-to-coarse" model that starts from the description of 
the changes of prices in continuous time and that allows one, from the mi- 
croscopic properties of the model to recover a large scale diffusion behaviour. 
More precisely, be defining a (multivariate) tick-by-tick model by means of 
marked point processes with appropriate stochastic intensities, we are able 
to control its features at all scale (i.e. its aggregation properties). Our 
model relies on multivariate Hawkes processes \n\ [TU] associated with posi- 
tive and negative jumps of the asset prices. Notably, by coupling suitably the 
stochastic intensities of upward and downward changes of prices for several 
assets simultaneously, we can reproduce microstructure noise {i.e. strong 
microscopic mean reversion) and the Epps effect (i.e. the decorrelation of 
the increments in microscopic scales) while preserving a standard Brownian 
diffusion behaviour on large scales. 

The paper is organized as follows: in Section [2] we define a univariate 
version of the model and compute the expected signature plot. We discuss 
how to simulate the model in practice and how to estimate its parameters 
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from real data. The vector (bivariate) version of the model is defined in 
Section [3l We show that a closed form expression for the price variations 
correlation functions can be explicitly obtained for all time scales and time 
lags. From such quantity one can deduce all second order properties of the 
process like individual signature plots, the Epps effect or lead-lag effects. 
We discuss the large scale diffusive limit of our model as far as its correla- 
tions properties are concerned and make the link with further mathematical 
asymptotic results (that are presented in a separate forthcoming paper). 
Comparisons to empirical data are provided in Section [S] in both ID and 2D 
frameworks. Section [6] is devoted to a conclusion and prospects for future 
research. 

2 The model in the univariate case 
2.1 Construction of the model 

We start with two point processes Ni{t) and N2{t) for t £ [0,T] that repre- 
sent respectively the sum of positive and negative jumps of some asset price 
X{t) over some time horizon [0, T]: 

X{t) = Ni{t) - N2{t). 

If Ni(t) and N2{t) are two independent Poisson processes with intensity 
it is easy to show that the model diffuses at large scale, i.e., when T — )• 
oo, by introducing the scaling factor T^/^, we obtain the following limit in 
distribution: 

lim -^X(tT) \/2ilB(t), t G [0,1], 

T-^+oo y/T V r- y n 

where \/2jI is the diffusive or macroscopic volatility, that accounts for the 
activity of negative and positive jumps, hence the factor 2 in the limit. 
According to Eq. ([T]), the corresponding mean signature plot is flat: for all 
r > 0, we have E[C(t)] = 2/i. In order to account for the previously reported 
noise microstructure features, intuitively and as confirmed by empirical ob- 
servations, one has to introduce some mean reversion in the small scales, 
while ensuring that this mean reversion effect vanishes on large scale. This 
can be naturally done within the context of (multivariate) Hawkes process 
[ni do] as foUows. 

For technical reason, we extend in a first step the time horizon [0, T] over 
the whole real line M = (— oo,-|-oo). Let Xi{t) for i E M be the stochastic 
intensities of two counting processes Ni{t), i = 1,2, such that at time t: 

\{t) = lim A-iE [N,{t + A) - N,{t) \ j;] (3) 
A--i>0 

where J-t stands for the filtration generated by the history of the processes 
Ni{t), N2{t). The bivariate process {Ni{t) , N2{t)} is a linear Hawkes process 
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if Ni(t) and N2{t) have no common jumps and if there exist fom' nonnegative 
functions {^ij}ij=i,2 such that 



Xi{t) = fii+ ipiiit - u)dNi{u) + ipij{t-u)dNj{u). (4) 

The so-obtained process can be shown to be well defined and to admit a 
version with stationary increments under the stability condition 

all the eigenvalues of the matrix {||99jj||i}are < 1, (5) 

where \\^p\\i = J^if{t)dt, see [T71 [10]. Mean reversion can be translated 
by the fact that the more X{t) goes up, the greater the intensity X2{t) 
and conversely, the more X{t) goes down, the greater the intensity Ai(t). 
This leads to the following simplified version of previous model (where only 
mean-reverting terms were kept): 

Ai(t) = fi+ f ip{t - s)dN2{s) (6) 

J — oo 

A2(t) = t Ht-s)dNi{s) (7) 



where is an exogenous intensity and ip{t) a positive kernel which is causal 
(i.e., Supp((^) C M"*"). Eqs. and ([7]) define two mutually exciting point 
processes that are stationary and stable under the condition ||v3||i < 1. 
A simple and natural choice for (p is a right-sided exponential function: 

ip{t) = ae-^H^+{t) (8) 

where a, /3 > are such that 

(9) 



/3 

2.2 Signature plot 

We are interested only in X(t) on [0, T] C M_|_ and for simplicity and without 
loss of generality, we se 10 X{0) = 0. Using Eq. i^, the mean signature plot 
can be written as 

C(r) = E[C(t)] = ^E[\X{in + l)r) - X{nT)\'] = ^EiXirf]. 

If {p{t) is defined as in Eq. then one can actually obtain a closed form for 
the mean signature plot using the approach initiated in OUT]. In Appendix 
1, we prove the following result: 



^Indeed, since only the increments of the processes Ni{t) come into play, we may (and 
will) assume that Ni{t) = N2{t) = hence X(0) = 0. 
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Proposition 2.1. Under the stability condition ([9]), we have 



C(r) = A + (1 - K^) 



1 - e-T^ 




(10) 



where 



A 



2/i 



1 



and 7 = a + /3. 



1 - llv'lli 



l + ll^lli 



One sees in particular a cross-over from tlie microstructural variance 



In Fig. 1(b), we display a plot example of the function C(r). 

2.3 Numerical simulations and parameter estimation 

In this section, we focus on the numerical simulations and the parameters 
estimation of the simplified univariate model defined by Eqs ^ and ([7|) in 
the case the function ip is the right sided exponential function Thus, 
there are 3 parameters, namely : 6 = {f-i,a,/3). 

Simulation of this process on an interval [0, T] can be performed using 
the thinning algorithm described in |27j . It basically consists in simulating 
on [0, r] a standard Poisson process with an intensity M large enough such 
that it satisfies the following condition: 



(which of course can be checked only a posteriori). A thinning procedure is 
then applied to each jump of the obtained process from the first one to the 
last one allowing to either reject the point (with probability M— Ai(t) — A2(t)) 
or mark it as a jump of Ni(t) with probability Ai(t) or of N2{t) with prob- 
ability X2{t), where t is the time of the considered jump. 

A realization of the process X over T = 42 hours is represented in Fig. 
dl^a) with e = {fi = 0.016, a = 0.023, /3 = 0.11). Let us note that n,a and /3 
are all expressed in the the same unit, namely seconds"^. These particular 
values were chosen to match the estimated parameters on real data (see 
Section 15. 2p . 

Let {t\^^}o<i<Ni{-) (resp. {tP}o<i<N2{-)) the upward (resp. downward) 
jumps of the realization of X. The estimation of the parameters can be 
processed in very different ways depending on what is the focus of the model. 
On the one hand, if one is mainly interested in the ability of the model to 



Vo = E[C7(0)] = A = 2E(Ai) 



to the diffusive variance 



Voo = E[d{oo)] = Ak^ 



\i{t) < M and X2{t) < M, for all t e [0,r] 
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reproduce the mean signature plot, the parameters can be estimated using a 
best fit of the reahzed signature plot. The realized signature plot over [0, T] 
is defined as 

T j T 

C{T) = ]-Y.\X{{n + \)T)-X{nT)\\ (11) 



T 

n=0 



and TcgrcssioTi estimator ^i-cg then naturally given by 

^reg = Argmin,|a(r)-C(T)|', (12) 

where C(r) is defined by (fTO|) . 

On the other hand, if the goal of the model is not simply to reproduce 
the signature plot behaviour but to mimick the arrival times themselves, 
it is more natural to consider the Maximum Likelihood Estimator (MLE) 
instead. This is possible since there is a closed formula for the likelihood. Let 
N(t) be a multivariate point process with conditionnal intensity A(i) = Xg(t) 
depending on a parameter 9. If, for every 9, the law of the process N(t) 
is absolutely continuous w.r.t. to the law of a standard stationary Poisson 
process, then the statistical model generated by the continuous observation 
of N{t) over [0, T] has a log-likelihood function given by 

L{9)= [ logXe{t)dN{t)+ [ {1-Xe{t))dt, (13) 
Jo Jo 

see for instance |2H [TO] or [26] in the context of Hawkes processes and the 
references therein. Consequently the likelihood functions 9 i— ?■ Li{9) and 
9 I—)' L2{9) generated by the observation of Ni{t) and N2{t) over [0,T] are 
given by 



0<if'<Ari{T) 0<tf'><N2iT) 



-/3(t(i)-t(.2') 



0<tf' <Af2(T) 



and similarly 



0<tf^<N2{T) 0<t^p<N2{T) 



0<tf'><NiiT) 



and the MLE of 9 based on the observation of X{t) over [0, T] is thus given 

by 

^MLE = ATgming{L,{9) + L2{9)). (14) 
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Of course, for both estimators 9,-cg and Omle the minimization (Eqs (jlip 
and (jl3p ) must be performed under the constraints 



/X > 0, a > 0, /? > 0, and the stabihty condition : — < 1. 



(15) 



Let us note that whereas the regression estimator ^reg can be performed on 
uniformly sampled data, the ML estimator Omle needs to have access to the 
point process itself. In that sense, Oreg can be considered as a (multi-scale) 
low-frequency estimator. For that reason, when applied to real data, we 
expect it to be much more stable than Omle (see Section [5]). Moreover, the 
regression estimator has the advantage to be computationally faster than 
the MLE. 

The estimated parameters of the realization shown in Fig. [U^a) using 
MLE are /2 = 0.016, a = 0.024 and /3 = 0.11. They match the true parame- 
ter values (/i = 0.16, a = 0.024, /3 = 0.11). 

The estimation using the regression estimator appears to be almost as 
accurate as the MLE: we obtain Ji = 0.016, a = 0.021 and /? = 0.010. 
In Fig. mb) we display the estimated and theoretical signature plots (Eq. 
()10p ). One sees that the latter provides a good fit of the data. 

3 The model in the bivariate case 
3.1 Definition 

The bivariate model is a natural extension of the former construction. We 
start from two processes Xi{t) and X2{t) with t G M, each constructed as 
in Secti0n r2.ll (i.e. in dimension 1) and introduce a supplementary coupling 
on the intensities of the processes in order to create a dependence struc- 
ture. More precisely, we consider four mutually exciting point processes 
{A'j(t)}j=i^,,,^4 associated with the positive and negative variations of Xi{t) 
and X2it): 



The joint law of the processes Ni{t) is characterized by their intensities: 



Let us note that, as previously, the kernels ipij account for the mutual and 
cross excitations of positive/negative parts of the couple of assets. In the 
sequel, in accordance with the univariate case, we suppose that /.ii = fi2 and 
= M4 (but we do not assume that = ns). 



Xi{t) 



Ni{t)-N2{t) 
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10 20 30 40 

time (h) 




50 100 150 200 

T (S) 

Figure 1: Numerical simulation of a ID price model using mutually exciting 
Hawkes processes with parameters ^ = 0.16, a = 0.024, /? = 0.11 (Let us 
note that these values were chosen to match the estimated parameters on 
real data as shown in Section [5^2]) . (a) Simulated sample path X(€) (42 hours 
long) (b) Estimated signature plot C(r) with the corresponding expected 
analytical shape according to Eq. ()10p . r is expressed in seconds. 
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In order to account for mean reversion and cross coupling between the 
two assets, we do not consider all possible mutual and cross excitations. 
More precisely, we only want to consider Upward- Xi-Upward-X2 and Downward- 
Xi-Downward-X2 couplings. We ignore further possible couplings Upward- 
Xi-Downward-X2 and Downward-Xi -Upward- between Xi and X2. We 
thus choose the matrix $ = {v3ij}i<ij<4 of the following form 



ifi2 v'ls 

Lp-iA 
\ 9931 V334 / 



Let us consider the matrix T = {rjj}i<j.j<4 with entries Tij 
and A the vector of mean intensities: 



Ai = E 



<m 

dt 



It is shown in Hawkes [17] that: 

A = (Id - ] 

The solution to this equation reads: 

Ml r34— /Ji— ri3 ^3 

ri2 r34— ri2— r34+i— r3i ri3 

Ml r34-/Ji-ri3 M3 

ri2 r34— ri2— r34+i— r3i ri3 

— MS ri2+r3i /J1+/J3 
ri2 r34— ri2— r34+i— r3i ri3 

-^3 ri2+r3i fii+fis 
ri2 r34— ri2— r34+i— r3i ri3 



(16) 



(17) 



(18) 



(19) 



As for the univariate case, we consider functions tpij that are causal 
exponentials: 

ipij{x) = a^je-^'i^'lK^ix). (20) 



In that case one simply obtains: 



a,; 



(21) 



For some to £ find a scale r, let us define the increments of a process 
X as 

ArX{to) = X{to + T)-X{to). 

If one wants to characterize multiscale self and cross-correlations of two 
asset prices Xi{t) and X2{t) the quantity of interest is the covariance matrix 
C(t,r) = {Cki{t,T)}i<k,i<2 with entries 



Ckiir, t) = Gov [ArXkito), A,X,(to + t)] . 



(22) 
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Thanks to the stationarity of the increments of X^, the matrix C does not 
depend on to- In particular, we recover from C(t,r) the mean signature 
plots on the diagonal and the Epps effect off diagonal as r varies. We also 
can estimate possible lead-lag effects across various time scales, that we 
can loosely define as the property that the function t i— t- Cki{T,t) is not 
symmetric around for k ^ I (see below). In the following, we denote by 
Cki{T) = Cki{T,t = 0) the correlations at lag t = 0. 

3.2 Computation of the signature plot and the Epps effect 

In Appendix 2, we show that we obtain a closed form expression for the 
covariance matrix (|22p . The general expression where all parameters are 
arbitrary is cumbersome and we limit ourselves to a matrix $ of the form 
([inD with (pij{x) = Oij exp{-l3ijx)lR^{x) and 

Pij = /3, i,j = 1,3. 
In that case we get the following proposition: 

Proposition 3.1. The expression of the covariance Ckiir) as a function of 
the time scale r reads: 

Cii{t) _ 2yliie-^i" 2Biie-'^^" ^ 2Aii 2Bii 2Aii 2Bii 
r Glr ^ GIt +^^^1+ Gi G\t GIt ' 

- - (^^^ + Glr + + + ^''^^ Glr + G^)' 

C22(t) _ 2^226"^!" 25226"^"" ^. 2^22 2B22 _ 2^ _ 2^22 

r Glr ^ Glr +^^^3+ + ^2^, 

with explicit constants provided in Appendix 2. 

It is interesting to discuss the fully symmetric case, i.e. when ipi2 = 9^34, 
9313 = (^31 and = = ^. In this case a direct computation from (jl9p 
leads to 

A = Ai = A2 = A3 = A4 = z „ (23) 

i — i 12 — i 13 

If follows from the results of Appendix 2, that 

-/i(r22 + ri2-ry 



Qi — Qa 

Q2 =Q3 



((ri2 + 1)2 - rl,) (1 - ri2 - ri3) 
-^ri3 



((ri2 + 1)2 - r?3) (1 - ri2 - ri3) 

After some algebra, one obtains: 

Corollary 3.2. In the fully symmetric case, the covariance matrix has the 
following expression: 

Cii(r) _ RCi RC2 
T 2Gi 2G2 

C2Gi^e"-G2 _ c^G2^ + QiG2^e-"gi - Cad^ 

^ 2G2^Gi^T 



12 



and 



Ci2(r) -RCi ^ RC2 



+ 



2Gi 2G2 



where 



R 

Ci 
C2 



2G2^Gx\ 

Pp 



ri2 + ri3-i 

(2 + ri2 + ri3)(ri2 + ri3) 

i + ri2 + ri3 
(2 + ri2-ri3)(ri2-ri3) 



i + r 



12 — i 13 



and 



Gi = /3(i + ri2 + ri3), G2 = /5(i + r 



12 — i 13 j 



One can reproduce the Epps effect by evaluating the behavior of the cor- 
relation coefficient p{t) = Ci2(t)/Cii(t) as a function of r. From Corollary 
()3.2p one can see that 

. , R{Q2 - Ql) , ^ /■ 2\ , n 

PiV = ^ T + 0{t ) as T ^ 0, 

R2 — Ri 

P(t) — )• — ; — as r — )• 00, 

' 2K + R1+R2 

where we have set Ri = RQ-jXi. If one considers the definition of each 
constant, one gets 

2ri3(i + ri2) 

'^'^ ^ i + rf3 + 2ri2 + rf, " ^ ^'^^ 

Notice that when the processes are not correlated, i.e. when ri3 = 0, we 
have Gi = G2 = /3(l+ri2) = /3+ai2 = 7, A = /^/(l-ris), R = /3^/(ri2-l), 
Q = = Q2 = ri2(ri2 + 2)/(l + ri2). Thus 

RQ _ ri2(ri2 + 2) 



Ai '^(ri2-i)(ri2 + i)2 



and therefore 



Gi G^T 



" ^ ^ +(1 ' 



12 



:i+ri2)2 ^ (l + ri2)2^ 7T 



We thus recover in that case the univariate result of Eq. 
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In order to account for the existence of a possible lead-lag effect, one can 
introduce a measure of the asymmetry of the covariance matrix at scale r 
as e.g. 

A(r) = Ci2(r,T)-C2i(r,T) . (25) 

This coefficient measures the difference of the correlation of assets 1 and 2 
returns at scale r and lag r when assets 1 is respectively in the past or in the 
future of asset 2. It is easy to see that if the matrix ()16p is symmetric then 
A = and there is no lead-lag. In the general case, from the results (and 
within notations) of Appendix 2, if Di = A12 — A21 and D2 = B12 — B21 
then: 

^ 1 + e-2G'i^2 - 2e-'^i^ ^ 1 + 6-2^2^2 _ 2e-<^2r 
^{r) = Di -2 + D2 -^2 (26) 

The expressions for Di and D2 are relatively heavy to handle in the general 
case. In the simple asymmetric case when 013 = one has: 

^ 2ai2Q3i/3m(2/3 + ai2) .^7) 

^ (ai2 - a34)(2/3 + 034 + ai2)(/3 - ai2) 

The specific study of lead-lags effects within this approach will be the scope 
of a forthcoming paper. 

3.3 Numerical simulations and parameter estimation 

In this section, we focus on numerical simulations and parameters estimation 
of the bivariate model in the fully symmetric case {'^12 = '■pj.i, ^12, = '■pj.i and 
A^i = = /^)i in the case where all the functions ipij are right sided expo- 
nential functions (|20p with = (3 for all i, j. Thus, there are 4 parameters, 
namely : 9 = (/i, ai,2, ai,3, /?)• 

A realization of this process bivariate {Xi,X2) over T = 20 hours is 
represented in Fig. [2] with Q12 = 0.23, an = 0.05, /3 = 0.11 and ^ = 0.015. 
(let us note that ^u, ai2, 013 and /3 are all expressed in the same unit, namely 
seconds"^). According to Eq. ()24p . the asymptotic correlation between the 
(large scale) increments of Xi and X2 is p ~ 0.15. Thus the two components 
are only moderately correlated (at a visual level, the graphs do not look 
alike). Fig. [3] shows the realization of the bivariate process (Xi,X2) in the 
case ai2 = 0.23, Q13 = 0.05, /3 = 0.11 and ^ = 0.015. Eq. ^ shows that 
the correlation between the two components is quite strong: p ~ 0.65. This 
can be clearly seen on the two graphs which look quite much alike. 

In Fig. m^b), the estimated and theoretical Epps effect (i.e. the function 
pij) = Ci2(t)/Cii(t)) are plotted for three sets of parameters, correspond- 
ing respectively to an asymptotic correlation coefficient of 0.15, 0.40 and 
0.65. 

As for the univariate case, one could perform maximum likelihood esti- 
mation of the parameters. Bivariate formula for the likelihood are easy to 
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obtained following the computations in Section 12.31 In the bivariate model, 
the quantities of interest are not only the signature plots of both Xi and 
X2 but also the Epps effect between Xi and X2. The regression estimation 
now takes the form 

dreg = Argmine(ai|Cii(r) - Cn{T)\^ + a2|C22(r) - C22(r)P 
+ ai2|Ci2(r)-Ci2(r))|2), 

where Ckiir) = Cki{T, 0) is defined as in ([22]) and 

1 

Cki{r) = -J2^rXk{nT)A^Xi{{n + l)T). (28) 

n=0 

The constants ai,a2 and ai2 are constant weights that are used to fix the 
relative minimization error of each term. Both estimators (^mle and ^reg) 
lead to quite accurate results (with the same magnitude of errors in both 
cases) . 




5 10 15 20 

time (h) 



Figure 2: Numerical simulation of a 2D symmetric price model using Hawkes 
processes with parameters ai2 = 0.23, ais = 0.01, /3 = 0.11 and fj, = 
0.015. (a) Sample path of Xi{t) (b) Sample path of X2{t). The asymptotic 
correlation between large scale increments is p = 0.15. 
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Figure 3: Numerical simulation of a 2D symmetric price model using Hawkes 
processes with parameters ai2 = 0.23, ai3 = 0.05, /3 = 0.11 and ^ = 0.015. 
(a) Sample path of Xi{t) (b) Sample path of X2{t). The processes appear 
more correlated then in Fig. [2] because the asymptotic correlation of large 
scale increments is p = 0.65. 

4 Diffusive (large scale) limit of the model 

It is interesting to discuss the large scale limit of the processes obtained by 
our construction. For that purpose, Let us define the following normalized 
processes for the general A^— variate version of the model: 

{xP{t)},=i...N = {^Xi{tT)},=i,„M for t G [0, 1]. (29) 

The question we would like to address concerns the existence and the prop- 
erties of the macroscopic limit A^-variate process {-^'^°°'' (t)}i=i...Ar obtained 
by considering the large scale limit T — )• oo. As far as the second order 
(correlation) properties of these processes are concerned, from Propositions 
12.11 and 13. H one sees that in both the univariate and the bivariate case, a 
diffusive limit exists. For instance, in the bivariate case, the limit process is 



16 







100 



200 
T 



300 400 



Figure 4: Epps effect for simulation samples. The estimated correlation 
coefficients are compared to the expected analytical curves (solid lines). 
From top top bottom the cross correlation between the two asset returns 
increases from 0.15 to 0.65. 



characterized by the following covariance matrix: 

2Au 2Bn A12+A21 I B12+B21 . 
C(2) = Gi G2 + G2 (30) 

' ' I B12+B21 2A22 _^ 2B22 ' 



Gi ^ G2 Gi ^ G2 

By using limit theorem for semimartingales, it is actually possible to rig- 
orously obtain a limiting process which is a multivariate Brownian motion 
with the appropriate covariance matrix. We describe and prove this results 
in details in a forthcoming paper. 



5 Comparison to empirical data 

5.1 Presentation of the data 

Let us recall that our goal is not to reproduce perfectly the features observed 
from empirical data since on the one hand, ultra high-frequency tick-by-tick 
series can be defined in various way {e.g. mid prices, transaction prices, ask 
or bid prices and so on) , each of them having slightly different properties and 
representing an arbitrary projection of the complex order book dynamics. 
On the other hand, it is well know that market intraday fiuctuations are not 
stationary and are characterized by strong long term heterogeneities. For 
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computation and parsimony purposes, our model is somehow naive since the 
exponential shape of the kernels (pij are somewhat arbitrary and the chosen 
symmetries are arbitrary and not necessarily suitable in practice. We rather 
aim at providing a first-brick model that can be considered as a tool to 
understand how the discrete nature of price variations at fine scales can be 
aggregated at coarser scales and lead to cross correlated Brownian diffusion 
processes. In that respect, comparisons to real data have to be interpreted 
rather at a qualitative level. 

The data that have been used in this paper consist in tick-by-tick last 
traded price time series, with, for each trade, the corresponding volume 
and a flag indicating whether the trade corresponds to a sell order or a 
buy order. The prices are either Eurex Euro-Bund futures contracts or 
Eurex Euro-Bobl futures contracts which correspond respectively to long- 
term (8.5-10.5 years) or medium-term (4.5 to 5.5 years) debt instrument 
issued by the Federal Republic of Germany. Euro-Bund and Euro-Boble 
are well known to have highly correlated price variations. So they are good 
candidates to study the Epps effect. Both open from (local-time) Sam to 
10pm, but we shall restrict to the most liquid period : 8am to 5:15pm. 
However, the liquidity (and the volatility) is highly seasonal during the day. 
Our model does not account for such a seasonality, so we shall restrict the 
data to intraday periods for which the underlying stationarity assumption 
is reasonable. Based on rough empirical considerations, we select the time 
period 9am to 11am. Moreover all computations have been made on last 
traded prices of buy orders only. Choosing sell orders would not have change 
the results, however, taking into account in the same time-series both buy 
and sell orders would lead to a highly bouncing artefact that shall not be 
able to be captured by our modelling approach. 

To summarize, the computations are made on two datasets : 

• Dataset I : 9am to 11pm from 11/01/2009 to 12/15/2009 (21 days) 
on the contract maturity 12/2009 

• Dataset II : 9am to 11pm from 06/01/2009 to 08/01/2009 (41 days) 
on the contract maturities 06/2009 and 09/2009 (for each day the most 
liquid maturity is selected). 

5.2 Signature plots and Epps effect 

Figure [5][a) shows the logarithm of the last traded price (only buy orders) 
of the Euro-Bund during a whole day 11/03/2009 (the contract maturity 
is 12/2009). The signature plot displayed in Fig. ^ih) has been computed 
using the dataset I described at the end of Section 15.11 Each day is con- 
sidered to be an independent and identically distributed realization of the 
same process. In solid line we have superimposed the signature plot obtained 
when we fit the mean theoretical signature plot by a regression model in the 
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univariate version of the model. The estimates we obtain are fi = 0.016, 
a = 0.024 and /3 = 0.11. We see that the curve associated with the model 
fits the data relatively well. Notice that the MLE method also provides fairly 
good results (in that case, we obtain /i = 0.014, a = 0.030 and f3 = 0.08). 
However, as explained previously, one expects the MLE to be less stable 
with respect to "noise" than the regression method on the mean signature 
plot and we empirically observe that its performance are worse than for the 
regression estimator. On other instances, they could lead to dramatically 
bad results and this is the reason why we discard the MLE estimator further 
on. 

Figure [6] shows the logarithm of the last traded price (only buy orders) of 
the Euro-Bund top) and Euro-Bobl (bottom) during a whole day 11/03/2009 
(maturity 12/2009). One can directly observe the return correlations be- 
tween the two assets. The large scale correlation coefficient we find is close 
to p = 0.77. In Figures m a) and[71^b) we plot the signature plots associated 
with the two assets while^c) displays the estimated Epps effect as measured 
by the correlation coefficient Ci2(t) at different scales. The computations 
were made using the dataset II described at the end of Section 15.11 The 
solid lines represent the fits according to the regression method of the bi- 
variate model. One can see that although significant discrepencies between 
empirical and fitted data are observed, given the simplicity of the model, 
one can consider that it captures fairly well both variance and covariance 
features of assets from small to large time scales simultaneously. 



6 Conclusion and prospects 

We have proposed in this paper a simple bivariate tick-by-tick price model 
based on Hawkes (self and mutually exciting) point processes. We have 
shown that closed form expressions can be obtained for its second order 
properties at all time scales. This allows one to recover major high frequency 
stylized facts, namely the signature plot behavior and the Epps effect. We 
have shown that the model is easy to simulate and can be estimated with 
a MLE method or using a moment method. When compared to real data, 
we have seen that the 2D model hardly accounts for the exact behavior 
of signature plots and correlation functions. However, as emphasized in 
section [5l the approach introduced in this paper has to be considered as a 
simple framework that allows one to address issues related to the relationship 
between fine and coarse scale properties of market dynamics. It can also 
be used as a simple tool in order to investigate intraday market features 
using few parameters that are easy to interpret. For instance, according 
to the univariate version of the model, the variance associated with the 
microstructure is 

2fx 



A 



1 
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Figure 5: (a) Last traded (only buy market orders) prices path of the Euro- 
Bund contract on the 3rd of November 2009 (contract maturity 12/2009) 
from Sam to 10pm. (b) Associated mean daily signature plot (computed on 
the dataset I described at the end of Section 15. ip . The solid line represents 
the fit using the regression estimator Oreg for the ID Hawkes model (Section 
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Figure 6: (a) Last traded (only buy orders) prices path of the Euro-Bobl 
contract on the 3rd of November 2009 (contract maturity 12/2009) from 
Sam to 10pm. (b) Same as (a) but for the Euro-Bund. 
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Figure 7: Self and cross correlation between Euro-Bund et Euro-Bobl re- 
turns as a function of the time scale r. The computations have been made 
using the on the dataset II described at the end of Section ISTTl (a) Signature 
plot associated with the Euro-Bund (b) Same as in (a) for the Euro-Bobl. 
(c) Epps effect between Euro-Bund and Euro-Bobl as measured by the co- 
variance between the two asset returns at different scales r. In (a), (b) 
and (c) the heavy lines represent a fit of the empirical curves according to 
the 2D Hawkes model (Proposition 13 . 1 ( ) . 
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Figure 8: (a) Histogram of values of x (as defined in Eq. (j3ip ) for the 
Bund front month contract from 5/14/2009 to 12/31/2009. (c) Same plot 
for the Bobl over the same period. The mean values for the two assets are 
respectively x = 0.29 and x = 0.36. 
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while the large scale diffuse volatility is 



One sees that the mean reversion (as measured by < ||<^||i < 1) is softened 
at large scale thanks to the diffusion, however, the influence of ip does not 
completely disappear. From the above tow expressions of the microscopic 
and the macroscopic variance, the large scale effect of mean reversion can 
be quantified by the function: 

which has a minimum around x = 1/3. In Fig. [8] we plot the histograms of 
estimated values of \\'-p\\i during a period of 6 months for Bund and Bobl 
front month contracts (from 05/14/2009 to 12/31/2009). It is striking to 
observe that both distributions are peaked (with rather large deviations how- 
ever) around x = 1/3. This suggests that the market chooses the strength of 
microstructural mean reversion in order to minimize its long term volatility. 

In a future work, we will address the questions related to lead-lag effects 
that can be easily accounted, for each time scale, within our model. It will 
also be interesting to consider non parametric estimates of the kernel shapes 
{(pij) along the line of the work of [HO] and to consider various questions 
related to the well known non stationarity and long-range correlations of 
the volatility. Let us stress that it is also tempting to bridge the parametric 
approach advocated in this paper to the recent work of Joulin et al. [22] 
where the authors studied and quantified the effects of exogeneous news with 
respect to the endogeneous noise on the jumps of stock prices. 



Appendix 1: Signature plot in the univariate case 

In the univariate case the shape of the signature plot can be directly com- 
puted for an exponential kernel ip as given by Eq. ([8]). Indeed, if one defines 



A = E 

it directly results from Eqs ([6]) 

A 



E 



dN2 



dt \ ~ dt 
and dZD that 



1 



P — a 



(32) 



(33) 



From the definition of the signature plot and by symmetry between the 
processes Ni{t) and N2{t), we have, for r > 0: 



C(t) = - (E [nUt]] - E [Ni{T)N2iT)]) 



(34) 
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Let Mij{t) denotes the continuous part of the covariance oidNi[t) and dNj[t) 
and M{t) = Mii{t) - Muit). From the definition of Ni{t) and N2{t), we 
then have: 



E[iVi(r)2] 
E[Ni{t)N2{t)] 



At + A^r^ + 



Jo 



Mii{u — v)dudv 

^0 

Mi2{u — v)dudv 



and 



C{t) = 2 a + r 



^0 



Mi2{u — v)dudv 



(35) 



In order to compute Mij(T) and M(r), one can first estimate the conditional 
mean with respect to J-f+r and then perform unconditional averages. Using 
this trick, from the definition of \i{t) and A2(i), one finds: 



Mii(r) = / ip{T - u)Mi2{u) du 

J — oo 

Mi2{t) = kf{T)+ I ip{T - u)Mu{u)du 

J —oo 

and consequently, M(r) satisfies: 

M(r) = -(/7(r)A + / ip{T - u)M{u) du . 



(36) 



If one seeks for an exponential solution M (x) = ae to this equation, one 
finds: 



-(a + /3) 
Aa{a + 2/3) 
2(a + /5) 



and finally, from ()35p . 
C(r) = 2A 



+ 1 



1 _ g-(Q+/3)T' 



(1 + 1)2/ (a + /3)r 



(37) 



Appendix 2: Correlation function in the multivari- 
ate case 

Let us define M(r) as the (continuous part of the) covariance matrix of 
individual point densities and A the mean intensity vector: 



A,; 



E 



dNi 



dt 



Gov 



dNi{to) dNj{to + u) 



dt 



dt 
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It is straightforward to show that C can be obtained in terms of M as: 

t+T 

du I dv Kpy{v — u) , (38) 



Cp.{T,t) = 

where the matrix K is 

Ku{z) 
Kuiz) 
K2i{z) 

K22{z) 

with J defined as: 



(Ai+A2)<5(z) + Jn(z) 
J2i{z) 

{A3 + A4)S{z) + J22{z) 



(39) 



Jii(z) = Mii(z) + M22iz) - Muiz) - M2i{z) 

Jl2{z)= Mi3 (z) + M24 (2) - M23 (z) - Mi4 (z) 

J2i{z) = Msiiz) + Mi2{z) - Mii{z) - Ms2{z) 

J22{z) = Mssiz) + Mu{z) - Msiiz) - M43 



(40) 



Along the same hne as in [T7], one can show that the matrix M satisfies 
the following integral equation: 



M{u) = ^{u) diag (A) + / ^{u- v)M{v) dv for u > 0. 



(41) 



The values for u < can be obtained thanks to the equality Mij{—v) = 
Mji(v) that is a direct consequence of the definition of M. 

Let us denote by F{s) the unilateral Laplace transform of F{u). By 
taking the Laplace transform of Eq. (|4ip . one gets the following set of 
coupled linear equation of M(s): 



Mij{s) = Aj^ijis) + ^ik{s) (Mkj{s) + MjkWik)) ■ 



(42) 



By using Eqs. ()40|) . it follows that the Laplace transform of Jp^ satisfies the 
following linear system: 

Ju{s) = -<fi2{s)Jii{s) + ipi3{s)J2iis) - (fuis) Ai + A2 + Jn (/3i2 ) 

+(^13(s)Jl2(/3l3) 

^12(5) = -(^12(5)^12(5) + (^13(5)^22(5) + (^13(5) A3 + A4 + J22(/3l3) 
-(^12(5)J2l(/3l2) 

^21(5) = -(^34(5)^21(5) + (^31 (5)Jll(5) + (^3l(s) Ai + A2 + Jll(/33l) 
-(^34(5)Ji2(/334) 

^22(5) = -(^34(5)^22(5) + (^31(5)^12(5) - (^34(5) A3 + A4 + J22(/334) 
+(^3l(s)J2l(/33i). 
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The solution to this system reads: 

/ Vl+Vl 11034+^13 '"3 

/ 1 



As) 



+'Pl2+iP34,+'P3i fW — tPSl 'Pis 

V2+V2 y34+'|gl3 V4 

l+</'12+'i334+¥'34 ¥'12— ¥'31 PlS 

¥'31 VI+V3+V3 ¥'12 

l+¥'12+¥'34+¥'34 ¥'12— ¥'31 ¥'13 

. ¥'31 V2+V4,+V4 ifiU 1 

\ l+¥'12+¥'34+¥'34 ¥'12-¥'31 V13 / 



(43) 



where we have denoted by v the vector 

/ -'Puis) (Ai + A2 + Jll(/3l2)) + ^13(s)^12(/3l3^ ^ 
(^13(5) ^As + A4 + J22(/3l3)) - (^12(5) J2l(/3l2) 
(^31 (s) (Ai + A2 + Jll(^3l)) - <^34{s)Ji2{P3a) 
y -(^34(5) (A3 + A4 + J22(/?34)) + '^3l(s)^2l(/53l) J 

This vector can be computed if by evaluating J(s) at s = /3i2, 
one now considers that the Laplace transform of <^ij{x) is 



(44) 

,/334. If 

(45) 



it is possible to compute the expression of each component Jap and, by 
the inverse Laplace transform, to obtain the correlation matrix Cap{t,T). 
Accounting for the symmetries in Eq. (jl6p . in Eq. (j44p . we now have 
4 constants to determine: Qi = Jii(/3), Q2 = </i2(/3), Q3 = </2i(/5) and 
Qa = <^22(/3)- If we solve the linear system obtained by evaluating ([13|) in 
s = /?, we obtain. 



Qi = 
Q2 = 



Ai (af3iai3(l + 034 +13 + 012) - 3/3ct34C«i2 - Pal^ - <i34'»12 - 2ai2/3^ - al^^Zi) + Pf^i^'lz 



"130:31(23 H 


- "12 


+ "34) - 33^012 - "12034 - 


"34"12 


- 33^034 


- 233 - 


43034012 — 


3a§^ - o?23 


A3(2/32c.i3 - 


- "31 


"13 +2ai23oi3 +"34012"13 


+ 0343 


"13) + Ai 


(oi3"ii 


— "123031 - 


- "31"34"12) 


"13031(23 H 


- "12 


+ "34) — 33^012 — "12O34 — 


"34"12 


- 33^034 


- 233 - 


43034012 — 


3"i4 -"223 


Ai(2/32c,3i - 


- "31 


"13 +20343031 +"340i2"31 


+ "123 


"31) + A3 


(oi3"31 


— 0343013 - 


- "13"34"12) 


"13031(23 H 


- "12 


+ "34) - 33^012 - "?2034 - 


"§4"12 


- 33^034 


- 233 - 


43034012 — 


- o?2'5 


A3 ("31" 


13(1 


f "12 + 3 + "34) — 33o34"i2 


-3"i4 


- 012034 


- 2Q343 


= - 0§40i2^ 


+ 3Al"i4 


"13031(23 H 


- "12 


+ "34) - 33^0:12 - "^2034 - 


"i4"12 


-33 = 034 


- 233 - 


43034012 — 


3"l4 - "223 



where the expressions for the Aj are provided in Eq. (|19p . One can now 
inverse the Laplace transforms in Eq. ()43p and and, from Eq. ()39p the 
expressions of the functions Kpy{t). This leads to 



Kii{t) 

Kl2{t) 
K2l{t) 
K22{t) 



2A35(t) + A22e-^i* + 5226"^^* 
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where 



Gi = Y + Z 
G2 = Y-Z 

Y = /? + -(ai2 + 034) 
Z = (ai2 - 034)^ + 4013031 
and the expressions of constants Aij, Bij read: 

^11 = - (Qi (2013031 + + 2^ai2 - "34012) 
-<32(ai2ai3 + 2Zai3 + "34013) 
+Ai (4013031 + 2o?2 + 4Zai2 - 2034012)) /(4Z) 

-Bii = - (<5i (034012 - 2013031 - 0^2 + 2Zai2) 
+^2(013012 - 2Zoi3 + 034013) + 

+Ai(20340i2 +4Zoi2 - 20i2 - 4oi303i))/(4Z) 
^12 = - (Q3(20i303i + 2Z012 + O12 - O34O12) 

-(54(013012 + 2Zoi3 + 034013) 

-A3(4Zoi3 + 2O34O13 + 20i20i3))/(4Z) 
-B12 = -(Q3(0340l2 - 2O13O31 + 2Z012 - O12) 

+(54(034013 + 012013 - 2Z013) 

+A3(20i30i2 - 4Zoi3 + 20i3034))/(4Z) 
^21 = ((5l(0l203i + 2Z031 + O34O31) 

+(52(012034 - 0I4 - 2O13O31 - 2Z034) 

+Ai (2012031 + 2034031 + 4Zo3i))/(4Z) 

-B21 = {Q2{c4a - O12O34 + 2O13O31 - 2Z034) 
-(5l(0i203i - O34O31 + 2Z031) 
+Ai(4Z03i - 2O12O31 - 203403i))/(4Z) 

^22 = ((54(012034 - 2O13O31 - 2Z034 - 0I4) 
+(53(031034 + 2Z031 + O31O12) 

+A3(2oi2034 - 20I4 - 4013031 - 4Z034))/(4Z) 

-B22 = ((54(2013031 - O12O34 + 0I4 - 2Z034) 
+(53(2^031 - O12O31 - O34O31) 

+A3(4oi303i - 4Z034 + 20I4 - 20i2 034))/(4Z). 
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By performing the double integral (|38p. one finally obtains: 

Cii(r) ^ 2yliie-°i" 25iie-°=" ^ 2Aa 2Bn _ 2^ _ 2Bn 
r Gfr Gir + ^^^^ + Gi G2 G?r Gir 

G2i(r) Gi2(r) _ , , ,, e'^^"-l , , , ^ ,, e-^="-l , 

?2 



- (A21 +Ai2)(— + p^) + (521+Sl2)(— +„ 



G22(r) ^ 2^226-'"^" 2^226-'-=" ^ 2^22 2^22 27^22 2^22 

r Gfr Gir +^^^3+ Gi G2 G?r G^r ' 
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